Impact of the electron-electron correlation on phonon dispersions: 
failure of LDA and GGA functionals in graphene and graphite 
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We compute electron-phonon coupling (EPC) of selected phonon modes in graphene and graphite 
using various ab-initio methods. The inclusion of non-local exchange-correlation effects within the 
GW approach strongly renormalizes the square EPC of the A'l K mode by almost 80% with respect 
to density functional theory in the LDA and GGA approximations. Within GW, the phonon slope 
of the A'l K mode is almost two times larger than in GGA and LDA, in agreement with phonon 
dispersions from inelastic x-ray scattering and Raman spectroscopy. The hybrid B3LYP functional 
overestimates the EPC at K by about 30%. Within the Hartree-Fock approximation, the graphene 
structure displays an instability under a distortion following the A'^ phonon at K. 

PACS numbers: 71.15.Mb, 63.20.kd, 78.30.Na, 81.05.Uw 
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The electron-phonon coupling (EPC) is one of the 
fundamental quantities in condensed matter. It deter- 
mines phonon-dispersions and Kohn anomahes, phonon- 
mediated superconductivity, electrical resistivity, Jahn- 
Teller distortions etc. Nowadays, density functional the- 
ory within local and semi-local approximations (DFT) 
is considered the "standard model" to compute ab-initio 
the electron-phonon interaction and phonon dispersions^. 
Thus, a failure of DFT would have major consequences in 
a broad context. In GGA and LDA approximations^, the 
electron exchange-correlation energy is a local functional 
of the charge density and the long-range character of the 
electron-electron interaction is neglected. These effects 
are taken into account by Green-function approaches 
based on the screened electron-electron interaction W, 
such as the GW method^. GW is considered the most 
precise ab-initio approach to determine electronic bands 
but, so far, it has never been used to compute EPCs nor 
phonon dispersions. The semi-empirical B3LYP func- 
tional^ partially includes long-range Hartree-Fock ex- 
change. B3LYP has been used to compute phonon fre- 
quencies but, so far, not the electron-phonon coupling. 

The electron-phonon coupling is a key quantity for 
graphene, graphite and carbon nanotubes. It deter- 
mines the Raman spectrum, which is the most com- 
mon characterization technique for graphene and nan- 
otubes^i^, and the high-bias electron transport in nan- 
otubes^. Graphene and graphite are quite unique sys- 
tems in which the actual value of the EPC for some 
phonons can be obtained almost directly from measure- 
ments. In particular, the square of the EPC of the highest 
optical phonon branch (HOB) at the symmetry K-point 
is proportional to the HOB slope near K^. The HOB 
K slope can be measured by inelastic x-ray scattering 
(IXS)^'^ or by the dispersion of the D and 2D lines as a 
function of the excitation energy in a Raman experiment 
_5 ^io, 11^12, 13 ^ A careful look at the most recent data sug- 



gests that the experimental phonon slopes (and thus the the phonon dispersion ud vs. q . The phonon dis- 



EPC) are underestimated by DFT^. The abihty of DFT 
(LDA and GGA) in describing the EPC of graphene was 
also questioned by a recent theoretical work^'*. 

Here, we show that: i) the GW approach, which pro- 
vides the most accurate ab-initio treatment of electron- 
correlation, can be used to compute the electron-phonon 
interaction and the phonon dispersion; ii) in graphite and 
graphene, DFT (LDA and GGA) underestimates by a 
factor 2 the slope of the phonon dispersion of the highest 
optical branch at the zone-boundary and the square of its 
electron-phonon coupling by almost 80%; iii) GW repro- 
duces both the experimental phonon dispersion near K, 
the value of the EPC and the electronic band dispersion; 
iv) the B3LYP hybrid functional^ gives phonons close to 
GW but overestimates the EPC at K by about 30 %; v) 
within Hartree-Fock the graphite structure is unstable. 

In Fig. [1] we show the phonon dispersion of graphite 
computed with DFTgga In spite of the general good 
agreement with IXS data, the situation is not clear for the 
HOB near K. In fact, despite the scattering among ex- 
perimental data, the theoretical HOB is always higher in 
energy with respect to measurements and the theoretical 
phonon slope (for the HOB near K) is underestimating 
the measured one. It is also remarkable that while the 
DFT K frequency is ~ 1300 cm^^, the highest measured 
is much lower at ^ 1200 cm^^. 

The dispersion of the HOB near K can also be obtained 
by Raman measurements of the graphene and graphite 
D-line (^1350 cm"-'^)^^. The D-line frequency lud de- 
pends on the energy of the exciting laser e^. Accord- 
ing to the double-resonance mode l^^i^^ , el activates a 
phonon of the HOB with momentum q=K-|-Aq along 
the K-M linc^ and energy fia;^). Aq is determined by 

eK-Aq,7r* - eK-Aq,7r cl - hujo/"^, whcrc ek,7r/ir' is 
the energy of the tt/tt* electronic state with momentum 
k. Thus, by measuring vs. cl and considering the 
electronic tt bands dispersion from DFT one can obtain 
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FIG. 1: Upper panel: Phonon dispersion of graphite. Lines 
are DFT calculations, dots and triangles are IXS measure- 
ments from Refs. [yla, respectively. Lower panel: phonon 
dispersion of graphene from DFT calculations. Dashed lines 
are obtained by subtracting from the dynamical matrix the 
phonon self-energy between the vr bands (aJ^ in the text). 



persion thus obtained is very similar to the one from 
IXS data and its slope is clearly underestimated by DFT 
(Fig. [21 upper panel). The same conclusion is reached 
by comparing the D-line dispersion vs. (directly 
obtained from measurements) with calculations (Fig. [51 
lower panel). Note that the dispersions of the Raman 
2D-linc^ is consistent with the dispersion of the D and 
thus in disagreement with DFT (LDA and GGA) as well. 

The steep slope of the HOB near K is due to the 
presence of a Kohn anomaly for this phonon''. In par- 
ticular, in Ref. 0, it was shown that the HOB slope is 
entirely determined by the contribution of the phonon 
self-energy between 7r-bands, Pq, to the dynamical ma- 



trix, I?q 



Wo 



•y/lJq/m is the phonon pulsation, where 



TO is the mass. For a given phonon with momentum q. 



P^ = ^j2^3-±5i^ (1) 

^Vfc ^ ek.TT — ek+q,7r* 



where the sum is performed on Nk wavevectors all over 
the Brillouin zone, £'(k+q)i,kj = (k + q, i|AVq|k, j) is the 
EPC, AVq is the derivative of the Kohn-Sham potential 
with respect to the phonon mode, |k, i) is the Bloch eigen- 
state with momentum k, band index i and energy ek,i. 
n{n*) identifies the occupied (empty) 7r-band. In Fig. [1] 
we show a fictitious phonon dispersion uJ^ obtained sub- 
tracting Pq from the dynamical matrix (uJ^ = ^ Bq^/m) 
for each phonon. The HOB is the branch which is mostly 
affected and, for the HOB, uJq becomes almost flat near 
K. Thus, DFT (LDA or GGA) fails in describing the 
HOB slope near K, slope which is determined by Pq. Pq 
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FIG. 2: (Color online). Upper panel: dispersion of the high- 
est optical phonon in graphite near K. Calculations are from 
DFT, or corrected to include GW renormalization of the 
electron-phonon coupling. Here, the DFT dispersion is ver- 
tically shifted by -40 cm~^ to fit measurements. Dots and 
triangles are IXS data from Refs. dj^, respectively. Squares, 
plus and diamonds are obtained from Raman data of Refs. 
llOlllljlla . respectively, using the double resonance mode l^'^'^^ . 
Lower panel: dispersion of the Raman D-line. 



is given by the the square EPC divided by 7r-band ener- 
gies. Thus, the DFT failure can be attributed to a poor 
description of the EPC or of the 7r-band dispersion. 

In graphene and graphite, it is known that standard 
DFT provides an underestimation of the tt and 7r*-band 
slopes of ~ 10 — 20 %^^'^^ . A very precise description 
of the bands, in better agreement with measurements, 
is obtained using GW— We thus computed the tt- 
bands with DFT (both LDA and GGA)if5 and GW^^ 
and compared with Hartree-Fock (HF)^ and BSLYP^^. 
Details are in |2ll The different methods provide band 
dispersions whose overall behavior can be described by a 
scaling of the tt energies^^. The different scaling factors 
can be obtained by comparing Ae^: the energy difference 
between the tt* and tt bands at the symmetry point M 
(L) for graphene (graphite). Ae^ is larger in GW than in 
DFT (Tab. [I|. Thus, inclusion of the GW correction to 
the electronic bands alone results in a larger denominator 
in Eq. [l] providing a smaller phonon slope and a worse 
agreement with experiments. The underestimation of the 
K phonon slope in DFT is, thus, due to the EPC. 

The EPC can be computed with linear response as, 
e.g., in Ref. but, at present, the use of this technique 
within GW is not feasible. Alternatively, the EPC as- 
sociated to a phonon mode can be determined by the 
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a) r-E„ 



b) K-A' 




FIG. 3: a, b): patterns of the r-E29 and K-A'i phonons of 
graphene. Dotted and dashed lines are the Wigner-Seitz ceUs 
of the unit-cell and of the \/3 x \/3 super-cell, c): Hartree- 
Fock equilibrium structure. 



variation of the electronic band energies by displacing the 
atoms according to the considered mode. In graphene, at 
K, there are doubly degenerate tt electronic states at the 
Fermi level. The HOB corresponds to the phonon at 
r and to the A[ at K. As an example, we consider the 
EPC associated to the r-E2g phonon and we displace 
the atoms according to its phonon pattern (see Fig. [3]). 
Following symmetry arguments^^, one can show that, in 
an arbitrary base of the two-dimensional space of the tt 
bands at K, the Hamiltonian is the 2x2 matrix: 



H = 2 



a b 

b* -a 



(2) 



1, and 



where each atom is displaced by d, |ap + |&p 

{D^)f — ^21^ I^Ki,Kj where the sum is performed 
on the two degenerate tt bands. Diagonalizing Eq. [21 
we see that an atomic displacement following the r-E2g 
phonon induces the splitting AEr — eK.ir* — ^k.tt and 



{Dl) 



1 I' AEr 

d^Q 16 V d 



^ ^ F = lim — 



(3) 



In analogous way, we define {Dj^)F = 
l-0(2K)i,KiP/4 for the A'^ phonon at K. Let 

us consider a ^/S x -s/S graphene supercell. Such a cell 
can be used to displace the atoms following the K-A'j^ 
phonon (Fig. [3]), since the K point is refolded in T. Let 
us call Ai?K the sphtting of the ck.tt bands induced 
by this displacement (since K is refolded in F, here 
ek.tt denotes the energies of the F band of the supercell 
corresponding to the tt band at K in the unit cell). 
Considering the atomic distortion of Fig. [31 displacing 
each atom by d, one can show that 



lim i 

d^o 8 \ d 



(4) 



In practice, by calculating band energies in the distorted 
structures of Fig. [3] and using Eqs. [3l [4] one obtains the 
EPCs of the r-E2g and K-A'j^ phonons between tt states. 
Similar equations can be used for graphite^"^. Results are 
in Tab.Utogether with the computed phonon frequencies. 
The EPCs from DFTqga are in agreement with those 
from linear response^. We also remark that, within the 
present " frozen-phonon" approach, the Coulomb vertex- 
corrections are implicitly included within GW. 



TABLE I: Electron-phonon coupling of the r-E25 and K-A'^ 
phonons computed with various approximations. Aeg (eV), 
{-Dq)p (eV^/A^) and Oq (eV/A^) are defined in the text, ujt 
(cjk) is the phonon frequency of the E23 (A'l) mode (cm~^). 
The GW wk for graphite (in parenthesis) is not computed 



directly (see the text). 



DFTlda 

DFTgga 

GW 

B3LYP 

HF 



DFTlda 
DFTgga 
GW 



4.03 
4.08 
4.89 
6.14 
12.1 



44.4 
45.4 
62.8 
82.3 
321 



-1 is the imaginary unit. 
Graphene: 



Qr 
11.0 
11.1 
12.8 
13.4 
26.6 



1568 
1583 

1588 
1705 



89.9 
92.0 
193 
256 
6020 



QK 

22.3 
22.5 
39.5 
41.7 

498 



1275 
1303 

1172 

960X1 



Graphite: 



Aeg 
4.06 
4.07 
4.57 



43.6 
44.9 
58.6 



10.7 
11.0 
12.8 



1568 
1581 



88.9 
91.5 
164.2 



21.8 1299 
22.5 1319 

35.9 (1192) 



To study the effect of the different computational 
methods on the the phonon slope (which is determined 
by Pq) we recall that Pq is the ratio of the square EPC 
and band energies (Eq. [1]). Thus, we have to compare 
ttq = {D'^)p / Aeg. As an example, assuming that the 
change of Pq from DFT to GW is constant for q near K, 



pGW 

q 

pDFT 

q 



„GW 



(5) 



and 



provides the change in the K phonon slope 



going from DFT to GW. To understand the results, we 
recall that in standard DFT the exchange-correlation de- 
pends only on the local electron-density. In contrast, 
the exchange-interaction in HF and GW is non-local. 
Furthermore, in GW, correlation effects are non-local 
since they are described through a dynamically screened 
Coulomb interaction. The hybrid functional B3LYP 
gives results intermediate between DFT and HF. 

Both ar and q:k are heavily overestimated by HF, the 
K-EPC being so huge that graphene is no more stable 
(the KA']^ phonon frequency is not real). Indeed, the 
HF equilibrium geometry is a 73 X V3 reconstruction 
with alternating double and single bonds of 1.40 and 
1.43 A lengths as in Fig. [31 (with a gain of 0.9 me V/ atom). 
These results demonstrate the major effect of the long- 
range character of the exchange for the K-EPC^'' but also 
the importance of the proper inclusion of the screening 
(included in GW but neglected in HF). Notice also that 
of graphite is smaller with respect to graphene by 



GW 



^--^10%. This is explained by the larger screening of the 
exchange in graphite (due to the presence of adjacent 
layers) than in graphene. On the contrary, within GGA 
and LDA, the graphite phonon frequencies and EPCs are 
very similar to those of graphene, since these functionals 
do not take into account the electron-electron interaction 
screening. 

Concerning the phonon slope, is 15% larger than 
Qfp^^. Indeed, DFT reproduces with this precision the 
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phonon frequency and dispersion of the HOB at T. On 
the contrary, is 60% larger than a^^^ , for graphite. 
This large increase with respect to DFT could explain the 
disagreement between DFT and the measured A'j^ phonon 
dispersion near K. To test this, we need to determine the 
GW phonon dispersion that, using Eq.[5]becomes ^ ~ 
^{BGW + rGWpDFTy^^ ^^^^^ ^GW ^ ^ q Moreover, 

we can assume ~ since the i?q component of 

the dynamical matrix (Eq. [T]) is not expected to have an 
important dependence on q (Fig.[T]). The value of B^^ 
is obtained as a fit to the measurements of Fig.[2F^. The 
resulting K A'j^ phonon frequency is 1192 cm^^ which 
is our best estimation and is almost 100 cm~^ smaller 
than in DFT. The phonon dispersion thus obtained and 
the corresponding D-line dispersion are both in better 
agreement with measurements (Fig. [5]). 

The partial inclusion of long-range exchange within 
the semiempirical B3LYP functional leads to a strong 
increase of the EPC at K as compared to the LDA and 
GGA functionals. However, comparing to the GW value, 
the EPC is overestimated by 30% and the corresponding 
frequency for the K.-A[ mode at 1172 cm~^ falls well be- 



low the degenerate K-mode which is around 1200 cm~^ 
in the experiment^ii^ (Fig- HI) and at 1228 cm~^ in our 
phonon calculation with B3LYP. We have checked that 
tuning the percentage of HF-exchange in the hybrid func- 
tional allows to match the EPC value of the GW ap- 
proach (in which case, the K-A'j^ mode remains the high- 
est mode. This may be a good way to calculate the full 
phonon dispersion of graphite/graphene within DFT, yet 
with an accuracy close to the one of the GW approach. 

Concluding, GW is a general approach to compute ac- 
curate electron-phonon coupling where DFT functionals 
fail. Such a failure in graphite/graphene is due to the 
interplay between the two-dimensional Dirac-like band 
structure and the long-range character of the Coulomb 
interaction-^''. However, GW can be also used in cases 
(in which the EPC is badly described by DFT) where 
the electron-correlation is short ranged2^. 
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In graphite, at the high-symmetry H point the four vr 
bands are degenerate two-by-two, Aeo being the energy 
difference. By displacing the atoms according to the T 
E2g phonon, these bands remain degenerate and the en- 



ergy difference is increased by Ae. In analogy to Eqs. |3]- 
Hwe define {D^)p = (Ae^ - Aeo)/(16d^). By displacing 
the atoms according to the K A^ phonon, the four bands 
are no longer degenerate, being n* (tt) the two bands 
which are up(down)-shifted and Ae — — e.^. We de- 
fine (D^)^ = (A? - Aeg)/(8d^), where A? indicates the 
average between the four possible tt* — vr couples. 
In principle, a direct calculation of cjk^ (and thus of 
Bk^) could be obtained, e.g., by finite difi'erences from 
a prohibitively expensive GW total energy calculation. 
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